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Abstract 

A Poisson-Fermi model is proposed for calculating activity coefficients of single ions 
in strong electrolyte solutions based on the experimental Born radii and hydra¬ 
tion shells of ions in aqueous solutions. The steric effect of water molecules and 
interstitial voids in the first and second hydration shells play an important role 
in our model. The screening and polarization effects of water are also included in 
the model that can thus describe spatial variations of dielectric permittivity, water 
density, void volume, and ionic concentration. The activity coefficients obtained by 
the Poisson-Fermi model with only one adjustable parameter are shown to agree 
with experimental data, which vary nonmonotonically with salt concentrations. 


1 Introduction 


Comprehensive discussions of theoretical and experimental studies on the ac¬ 
tivity coefficient of single ions in electrolyte solutions have been recently given 
by Fraenkel |T], Valiko and Boda [2], and Rowland et al. |3], where more refer¬ 
ences can also be found. The Poisson-Fermi (PF) model proposed in this paper 
belongs to the continuum approach that traces back to the simple, elegant, 
but very coarse theory — the Debye-Hiickel (DH) theory. As mentioned by 
Fraenkel, the continuum theory has evolved in the past century into a series of 
modihed Poisson-Boltzmann (PB) equations that can involve an overwhelm¬ 
ingly large number of parameters in order to £t Monte Carlo (MC), molecular 
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dynamics (MD), or experimental data. Many expressions of those parameters 
are rather long and tedious and do not have clear physical meaning [1]. 


The Debye-Hiickel model is derived from a linearized PB equation [1]. Ex¬ 
tended from the DH model, the Pitzer model |5] is the most eminent approach 
to modeling the thermodynamic properties of multicomponent electrolyte so¬ 
lutions due to its unmatched precision over wide ranges of temperature and 
pressure [5]. However, the combinatorial explosion of adjustable parameters 
in the extended DH modeling functions (including Pitzer) can cause profound 
difficulties in fitting experimental data and independent verihcation because 
the parameters are very sensitive to numerous related thermodynamic prop¬ 
erties in multicomponent systems [3] . The Poisson-Fermi model proposed here 
involves only one adjustable parameter. 


The ineffectiveness of previous Poisson-Boltzmann models is mainly due to 
inaccurate treatments of the steric and correlation effects of ions and water 
molecules whose nonuniform charges and sizes can have signihcant impact on 
the activities of all particles in an electrolyte system. Unfortunately, the point 
charge particles of PB theories have electric fields that are most approximate 
where they are largest, near the point. PB theories are not an appealing choice 
for the leading terms in a series of approximations, for that reason. The PF 
theory developed in our papers piTIIHII^IITU] demonstrates how these two effects 
can be described by a simple steric potential and a correlation length of ions. 
The parameters of the PF theory describe distinct physical properties of the 
system in a clear way [9]. The Gibbs-Fermi free energy of the PF model reduces 
to the classical Gibbs free energy of the PB model when the steric potential 
and correlation length are omitted [9]. The PF model has been verihed with 
either MG, MD or double layer data at (more or less) equilibrium piTIIS] . and 
nonequilibrium data from calcium and gramicidin channels inn]. 


Here, we apply the PF theory to study the activity properties of individual ions 
in strong electrolytes. The steric effect of all particles and the interstitial voids 
that accompany them are described by a Fermi-like distribution that dehnes 
the water densities in the hydration shell of a solvated ion and the particle 
concentrations in the solvent region outside the hydration shell. The resulting 
correlations produce a dielectric function that shows variations in permittivity 
around the solvated ion. The experimental concentration-dependent dielectric 
constant model proposed in [2] is used to dehne the concentration-dependent 
Born radii of the solvated ion in the present work. The experimental data of 
the activity coefficients of NaGl and GaGU electrolytes reported in m are 
used to test the PF model. 
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2 Theory 


The activity coefficient 7 * of an ion of species i in electrolyte solutions describes 
the deviation of the chemical potential of the ion from ideality ( 7 j = 1). 
The excess chemical potential is = ksTln.'-fi, where ks is the Boltzmann 
constant and T is an absolute temperature. In Poisson-Boltzmann theory, the 
excess chemical potential can be calculated by [12] 

= AGf® - AG“. AGf® = AG" = (1) 

where the center of the hydrated ion (also denoted by i) is set to the origin 
0 for convenience in the following discussion and qt is the ionic charge. The 
potential function of spatial variable r is found by solving the Poisson- 

Boltzmann equation 

-e,VV™(r) = ^ qjCj{v) = p(r), (2) 

i=i 

Cj{r) = Cf exp (-/lj0™(r)) , (3) 

where the concentration function Cj{r) is described by a Boltzmann distribu¬ 
tion (3) with a constant bulk concentration C^, = e^eo, is the dielec¬ 

tric constant of bulk water, and eo is the vacuum permittivity. The potential 
0°(r) of the ideal system is obtained by setting p(r) = 0 in (2), i.e., all ions 
of K species in the system do not electrostatically interact with each other 
since Qj = 0 for all j. We consider a large domain hi of the system in which 
0 PB(r) = 0 on the boundary of the domain dVl. The ideal potential ^^(r) is 
then a constant, i.e., AG° is a constant reference chemical potential indepen¬ 
dent of cf. 


For an equivalent binary system, the Debye-Hiickel theory simplifies the calcu¬ 
lation by analytically solving a linearized equation of (2) so that the potential 
function 0^®(r) becomes a constant [1] 


Qji^ f _ ( esksT \ 
dvre/ K \Tfj=iq]CfL) 


(4) 


dependent of the bulk concentration Cf, where L is the Avogadro constant. 


The Poisson-Fermi equation proposed in [5] is 

K+l 


j=i 

Cj(r) = Cf exp (-/?j0P^(r) + = In 


pB ’ 


(5) 

( 6 ) 
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where S'*’'^(r) is called the steric potential, r(r) = 1 — VjCj{v) is a void 

fraction fnnction, = 1 — VjCj is a constant void fraction, and Vj is 

the volume of a species j particle (hard sphere). Note that the PF equation 
includes water as the last species of particles with the zero charge qx+i = 0 . 
The polarization of the water and solution is an output of the theory. The 
water can be described more realistically, for example, as a quadrupole in 
later versions of the theory. The distribution ( 6 ) is of Fermi type since all 
concentration functions are bounded above, i.e., Cj{v) < 1/vj for all particle 
species with any arbitrary (or even infinite) potential 0 (r) at any location r 
in the domain hi [U]. The Boltzmann distribution (3) would however diverge 
if 0(r) tends to inhnity. This is a major deficiency of PB theory for modeling 
a system with strong local electric helds or interactions. The PF equation (5) 
and the Fermi distribution reduce to the PB equation (2) and the Boltzmann 
distribution (3), respectively, when Ic = = 0, i.e., when the correlation 

and steric effects are not considered. 

If the correlation length Ic = 2aj 7 ^ 0, the dielectric operator e = es(l — 
approximates the permittivity of the bulk solvent and the linear response of 
correlated ions [61I7II13II14] . where Oj is the radius of the ion. The dielectric 
function e(r) = 6^/(1 + rj/p) is a further approximation of e. It is found by 
transforming (5) into two second-order PDFs [B] 



(7) 

( 8 ) 


V^0^^(r) = ^/(r) 


by introducing a density like variable T that yields a polarization charge den¬ 
sity 7] = —CsT — p of water using Maxwell’s hrst equation [7]. Boundary 
conditions of the new variable T on the boundary dQ were derived from the 
global charge neutrality condition [ 6 ]. 

To obtain more accurate potentials at the origin 0, i.e., 0^^(O), we need to 
consider the size and hydration shell of the hydrated ion i. The domain hi is 
partitioned into three parts such that hi = hlionUnshUhlsoiv, where hlion is the 
spherical domain occupied by the ion i, hlsh is the hydration shell of the ion, 
and fisoiv is the rest of the solvent domain as shown in Fig. 1. The radii of hlion 
and the outer boundary of flsh are denoted by and respectively, 

whose values will be determined by experimental data. It is natural to choose 
the Born radius as the radius of flion [I2]- We consider both first and 

second shells of the ion [T5IIT6] . The dielectric constants in flion and r2\r2ion 
are denoted by eion and e*, respectively. 
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Fig. 1. The model domain is partitioned into the ion domain rijon (with radius 
shell domain flsh (with radius and solvent domain risoiv 

The PF equation (5) then becomes 


e - l) VV‘’''(r) = P(r) 


qi6{r - 0) in fijon 
QjCjir) in 


where 5(r — 0) is the delta function at the origin, /c = 0 in ffion, /c 7 ^ 0 in 
ff\nion, e = eionCo in f^ion, and e = = e^eo in f2\f2ion- The shell radius Rf^ 

is determined by Eq. (6) as 


ctrc 

^Sh 


, b'sh — 

In -— 

V^shFB 


In 


07 


^ — 


TB 

— O'^ 

cs * 


+ UwO' 


( 10 ) 


where is the volume of a water molecule and Vsh is the volume of the hydra¬ 
tion shell that depends on the bulk void fraction T®, the bulk water density 
C^, and the total number O 7 (coordination number) of water molecules oc¬ 
cupying the shell of the hydrated ion i. Note that the shell volume Vsh varies 
with bulk ionic concentrations Cf. The occupancy number O 7 is given by 
experimental data and so is the shell volume that of course determines 

the shell radius 


To deal with the singular problem of the delta function 5(r — 0) in Eq. (9), 
we use the numerical methods proposed in [6] to calculate 0^^(r) as follows: 

(i) Solve the Laplace equation V^0^(r) = 0 in f2ion with the boundary con¬ 
dition 0^(r) = 0*(r) = gi/(47reioneo k - 0|) on_5hlion. 

(ii) Solve the Poisson-Fermi equation (9) in with the jump condition 
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eV0^^(r) ■ n = 
condition 0^^(r) 
d^ion [ 6 ]. 


—eioneoV(0*(r) + 0'"(r))-n on c?f2ion and the zero boundary 
= 0 on dfl, where [m] denotes the jump function across 


The evaluation of the Green function 0*(r) on 5f2ion always yields hnite num¬ 
bers and thus avoids the singularity. Note that our model can be applied to 
electrolyte solutions at any temperature T having any arbitrary number (K) 
of ionic species with different size spheres and valences. 


3 Results 


Numerical values of model notations are given in Table 1, where the occupancy 
number O" = 18 is taken to be the experimental coordination number of the 
calcium ion Ca^^ given in [15] for all ions i = Na"*", Ca^^, and Cl“ since the 
electric potential produced by the solvated ion diminishes exponentially in the 
outer shell region in which a small variation of for i = Na"*" and Cl“ does not 
affect numerical approximations too much. Obviously the coordination number 
may be different for different types of ions and at different concentrations and 
so on. We were surprised that we can fit experimental data so well using a single 
experimentally determined occupancy number for all ions and conditions. 

As discussed in [2|, the solvation free energy of an ion i should vary with 
salt concentrations and can be expressed by a dielectric constant e{Cf) that 
depends on the bulk concentration of the ion Cf. Following |2], we assume 
that 

t(Cf) = £. - + {CfY'" (11) 

with only one parameter 6 i, whose value is given in Table 1, instead of two 
in [2|. Note that e{Cf) is a constant when the dimensionless Cf is given. It 
is not a function of a spatial variable r like e(r). The parameter hj represents 

the ratio of the factor of Cf to that of (^Cf^ in the original formula, where 
the factors of various electrolytes are taken from various sources of either 
theoretical or experimental data [2|. Our ratios Si in Table 1 are comparable 
with those given in [2|. 


The Born formula of the solvation energy can thus be modihed as 




Qi f± 
87reo0(Cf)Ji^ [e„ 





e(Cn (W - 1) 

^w(6(Gf)-l)’ 


( 12 ) 


where is the Born radius when Cf = 0 (6*(0) = 1) and = d(Cf)Rf 

is the concentration-dependent Born radius used to define Oion in Fig. 1 when 
cf 7 ^ 0. The Born radii in Table 1 are cited from [2], which are computed 
from the experimental hydration Helmholtz free energies of these ions given in 
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ini- All values in Table 1 are either physical or experimental data except that 
of 6 i, which is the only adjustable parameter in our model. All these values 
were kept hxed throughout calculations. 


Table 1. Values of Model Notations 


Symbol 

Meaning 

Value 

Unit 

ks 

Boltzmann constant 

1.38 X 10-23 

J/K 

T 

temperature 

298.15 

K 

e 

proton charge 

1.602 X 10-^® 

C 

eo 

permittivity of vacuum 

8.85 X 10-^^ 

F/cm 

^ion? 

dielectric constants 

1, 78.45 


Ic 

correlation length 

i = Na+,Ca2+, Cr 

A 

®Na+) ®Ca^+ 

radii 

0.95, 0.99 

A 

Oci-) ®H20 

radii 

1.81, 1.4 

A 

pO pO pO 

-^Na+’ -^Ca2+’ -^cr 

Born radii in Eq. (12) 

1.617, 1.706, 2.263 

A 

^Na+) ^Ca?+i ^Cl- 

in Eq. (11) 

4.2, 5.1, 3.8 


or 

in Eq. (10) 

18 



The PF results of Na+, Ca^’*', and Cl activity coefficients agree well with 
the experimental data P as shown in Figs. 2 and 3 for NaCl and CaCh 
electrolytes, respectively, with various [NaCl] and [CaCh] from 0 to 2.5 M. In 
Fig. 4, we observe that the Debye-Hiickel theory oversimplihes the Ca^^ activ¬ 
ity coefficient to a straight line as frequently mentioned in physical chemistry 
texts |1] because the theory does not account for the steric and correlation 
effects of ions and water, let alone the atomic structure of the ion and its 
hydration shell as shown in Fig. 1. Both PB and PF results in Fig. 4 were 
obtained using the same atomic Fermi formula (10) for shell radii Rf^ in fIsh 
and the same concentration-dependent Born formula (12) for Born radii 
in flion- Therefore, the only difference between PB and PF is in flsoiv, where 
Ic = = 0 for PB and C 7 ^ 0 and ^ 0 for PF. Note that these two 

formulas are not present in previous PB models. Fig. 4 shows that the corre¬ 
lation and steric effects still play a signihcant role in the solvent domain flsoiv 
although the domain is Rq’^ 2 + = 4.95 A (not shown) away from the center of 
the Ca^"*" ion. The ion and shell domains are the most crucial region to study 
ionic activities. For example, Fraenkel’s theory is entirely based on this region 
— the so-called smaller-ion shell region |T]. 

The PF model can provide more physical details near the solvated ion (Ca^’*', 
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O Experiment (Na"^) 



([NaCI]/M)^^^ 


Fig. 2. Comparison of PF results with experimental data m on i = Na^ and Cl 
activity coefficients 7 * in various [NaCl] from 0 to 2.5 M. 



Fig. 3. Comparison of PF results with experimental data [H] on i = Ca^"*" and Cl 
activity coefficients 7 * in various [CaCl 2 ] from 0 to 2 M. 

for example) in a strong electrolyte ([CaCl 2 ] = 2 M) such as the dielectric 
function e(r) of varying permittivity (shown in Fig. 5), variable water density 
Ch 2 o(i') (in Fig- 5), concentration of counterion (C'Q-(r) in Fig. 6), electric 
potential (0^^(r) in Fig. 6), and the steric potential (S'*’’^(r) in Fig. 6). Note 
that the dielectric function e(r) is an output, not an input of the model. 
The steric effect is small because the conhguration of particles (voids between 




0.5 



I_•_I_I_ 

0 0.5 1 1.5 

([CaCIJ/M)''^ 


Fig. 4. Comparison of Poisson-Fermi (PF), Poisson-Boltzmann (PB), and De- 
bye-Hiickel (DH) results on i = Ca^+ activity coefficients 7 * in various [CaCl 2 ] 
from 0 to 2 M. 

particles) does not vary too much from the solvated region to the bulk region. 
However, the variation of mean-held water densities C'H 2 o(r) has a signihcant 
effect on the dielectrics in the hydration region as shown by the dielectric 
function e(r). The strong electric potential in the Born cavity fljon and 

the water density C'H 2 o(r) in the hydration shell ffsh are the most important 
factors leading the PF results to match the experimental data. PF theory deals 
well with the much more concentrated solutions in ion channels where void 
effects are important [5]. 


4 Conclusion 


We have proposed a Poisson-Fermi model for studying activities of single ions 
in strong electrolyte solutions. The atomic structure of ionic cavity and hydra¬ 
tion shells of a solvated ion is modeled by the Born theory and Fermi distribu¬ 
tion using experimental data. The steric effect of ions and water of nonuniform 
sizes with interstitial voids and the correlation effect of ions are also consid¬ 
ered in the model. With only one adjustable parameter in the model, it is 
shown that the atomic structure and these two effects play a crucial role to 
match experimental activity coefficients that vary nonmonotonically with salt 
concentrations. 
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Fig. 5. Dielectric e(r) and water density C'H 2 o(r) profiles near the solvated ion Ca^"*" 
with [CaCl 2 ] = 2 M. 



Fig. 6. Cl concentration electric potential and steric potential 

gtrc(j,) pj-Qfiieg near the solvated ion Ca^"'' with [CaCl 2 ] = 2 M. 
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